
clear all
************************************
*Numeric Export disabled
local TEX ""
local TEXON 0 //Change "0" to "1" to enable exporting numeric results to .txt files
if !`TEXON' local TEX "*"
************************************
************************************
*Figure Export disabled
local FIG ""
local FIGON 0 //Change "0" to "1" to enable exporting figures (RELEVANT FOLDERS MUST FIRST BE CREATED)
if !`FIGON' local FIG "*"
************************************

**************************************************
**************************************************
**Evaluate the ability of analyses limited to preference gaps to 
**recover DGP under varios scenarios
**************************************************
**************************************************

**************************************************
*Must set working directory if not already set
**************************************************

******************************************************************
******************************************************************
**SAVE PARAMETERS GENERATED IN PROGRAM BELOW
**Results will post to: results_simngint.txt
******************************************************************
******************************************************************
use "Enns_SociologicalScience_Reproduction\DS1_0.dta", clear

//Recode outcome so dichotomous
gen policyadopted =.
recode policyadopted .=0 if OUTCOME==0
recode policyadopted .=1 if OUTCOME==2
recode policyadopted .=1 if OUTCOME==3
recode policyadopted .=1 if OUTCOME==4

//Recode variables as log of odds -- Gilens 2012, p.74
quietly: gen  pred90_sw_lor = log(pred90_sw/(1- pred90_sw))
label var pred90_sw_lor "log of odds: pred90_sw"

quietly: gen  pred10_sw_lor = log(pred10_sw/(1- pred10_sw))
label var pred10_sw_lor "log of odds: pred10_sw"

//Generate group intercept variable to account for the higher probability
//of policy adoption when the affluent support a policy more than other income
//groups (see, e.g., Figure 2 in the text)
**9010
quietly: gen dif_90_10=.
quietly: recode dif_90_10 .=1 if (pred90_sw - pred10_sw)>.1
quietly: recode dif_90_10 .=0 if (pred10_sw - pred90_sw)>.1

**********************
*Estimate parameters to be used in DGP for simulations below.
**********************
*90th income percentile and 10th income percentile
**********************

*90: Policy responsiveness to the 90th percentile
quietly: logit policyadopted pred90_sw_lor dif_90_10 if abs(pred90_sw-pred10_sw)>.1
//Save the coefficients from the previus model to be used in simulations.
mat betaDGP = e(b)

local coef_9010 = betaDGP[1,1]
local coefgap_9010 = betaDGP[1,2]
local constant_9010 = betaDGP[1,3]
*Divide coef by 2 because will be using it to generate responsiveness
*to both high and low income groups
local coefntdgp = `coef_9010'/2


*10: Policy responsiveness to the 10th percentile
quietly: logit policyadopted pred10_sw_lor dif_90_10 if abs(pred90_sw-pred10_sw)>.1
//Save the coefficients from the previous model to be used in simulations.
mat betaDGP = e(b)

local coef_109010 = betaDGP[1,1]
local coefgap_109010 = betaDGP[1,2]
local constant_109010 = betaDGP[1,3]
local coefntdgpt = `coef_109010'/2 //divde by 2 as explained above

*Save parameters used in dgp
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coefnt) result(`coef_9010') replace //beta (responsiveness to 90th)
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coefgapnt) result(`coefgap_9010') append //group intercept
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(consnt) result(`constant_9010') append //alpha constant
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coefntdgp) result(`coefntdgp') round(3) append //beta/2
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coeftnt) result(`coef_109010') append //beta (responsiveness to 10th)
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coefgaptnt) result(`coefgap_109010') append //group intercept
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(constnt) result(`constant_109010') append //alpha constant
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(coefntdgpt) result(`coefntdgpt') round(3) append //beta/2
******************************************************************
******************************************************************
******************************************************************
******************************************************************

******************************************************************
******************************************************************
**THIS PROGRAM IS THE SAME CODE AS ABOVE
**Can't save parameters for .tex file within program, so above is just to store parameters for .tex file
******************************************************************
******************************************************************
program sim_90_10

use "Enns_SociologicalScience_Reproduction\DS1_0.dta", clear

//Recode outcome so dichotomous
quietly: gen policyadopted =.
quietly: recode policyadopted .=0 if OUTCOME==0
quietly: recode policyadopted .=1 if OUTCOME==2
quietly: recode policyadopted .=1 if OUTCOME==3
quietly: recode policyadopted .=1 if OUTCOME==4

//Recode variables as log of odds -- Gilens 2012, p.74
quietly: gen  pred90_sw_lor = log(pred90_sw/(1- pred90_sw))
label var pred90_sw_lor "log of odds: pred90_sw"

quietly: gen  pred10_sw_lor = log(pred10_sw/(1- pred10_sw))
label var pred10_sw_lor "log of odds: pred10_sw"


//Generate group intercept variable to account for the higher probability
//of policy adoption when the affluent support a policy more than other income
//groups (see, e.g., Figure 2 in the text)
**9010
quietly: gen dif_90_10=.
quietly: recode dif_90_10 .=1 if (pred90_sw - pred10_sw)>.1
quietly: recode dif_90_10 .=0 if (pred10_sw - pred90_sw)>.1

**********************
*Estimate parameters to be used in DGP for simulations below.
**********************
*90th income percentile and 10th income percentile
**********************

*90: Policy responsiveness to the 90th percentile
quietly: logit policyadopted pred90_sw_lor dif_90_10 if abs(pred90_sw-pred10_sw)>.1
//Save the coefficients from the previus model to be used in simulations.
mat betaDGP = e(b)

local coef_9010 = betaDGP[1,1]
local coefgap_9010 = betaDGP[1,2]
local constant_9010 = betaDGP[1,3]
*Divide coef by 2 because will be using it to generate responsiveness
*to both high and low income groups
local coefntdgp = `coef_9010'/2


*10: Policy responsiveness to the 10th percentile
quietly: logit policyadopted pred10_sw_lor dif_90_10 if abs(pred90_sw-pred10_sw)>.1
//Save the coefficients from the previus model to be used in simulations.
mat betaDGP = e(b)

local coef_109010 = betaDGP[1,1]
local coefgap_109010 = betaDGP[1,2]
local constant_109010 = betaDGP[1,3]
local coefntdgpt = `coef_109010'/2 //divde by 2 as explained above

//Save the mean, min, and max of the policy prefernces of the 90th, 50th, and 10th income
//percentiles to be used in simulations 

quietly: sum pred90_sw 
local high_mean = `r(mean)'
local high_min = `r(min)'
local high_max = `r(max)'

quietly: sum pred10_sw 
local low_mean = `r(mean)'
local low_min = `r(min)'
local low_max = `r(max)'

quietly: sum pred50_sw 
local mid_mean = `r(mean)'
local mid_min = `r(min)'
local mid_max = `r(max)'

//Create a matrix with the mean values from above.
matrix mean = (`high_mean', `mid_mean', `low_mean')

//Generate coveriance matrix
quietly: cor pred90_sw pred50_sw pred10_sw, cov
mat cov = r(C)



*****************************************
//Generate simulated preferences of the 90th, 50th and 10th income percentiles
//based on characteristics of the observed data.
//Table 1 only analyzes 10th and 90th income percentiles, but 50th include
//to make sure data match observed characteristics of full analysis.

clear

*high represents 90th income percentile
*mid repersents 50th income percentile
*low represents 10th income percentile
drawnorm high mid low, n(2000) means(mean) cov(cov)

drop if high<`high_min' | high>`high_max' | low<`low_min' | low>`low_max'| mid<`mid_min' | mid>`mid_max'

//rescale preference variables.
quietly: gen logodds_high=ln(high/(1-high))
quietly: gen logodds_mid=ln(mid/(1-mid))
quietly: gen logodds_low=ln(low/(1-low))

//generate indicator for whether afflluent are more supportive or less supportive
//than low income.
quietly: gen dif_90_10=.
quietly: recode dif_90_10 .=1 if (high - low)>.1
quietly: recode dif_90_10 .=0 if (low - high)>.1


//generate predicted policy adoption based on actual parameters in the data above.

***************************************************************************************
*predictedoutcome9010_9010: responsiveness to both 90th & 10th income percentile
quietly: gen yhat_9010 = exp(`constant_109010' + (`coefntdgpt' * logodds_high) + (`coefntdgpt' * logodds_low) + (`coefgap_109010' * dif_90_10)) / (1 + exp(`constant_109010' + (`coefntdgpt' * logodds_high) + (`coefntdgpt' * logodds_low) + (`coefgap_109010' * dif_90_10))) if dif_90_10!=.
quietly:     gen predictedoutcome9010_9010 = 0 if dif_90_10!=.
quietly:     replace predictedoutcome9010_9010 = 1 if runiform() < yhat_9010 & dif_90_10!=.
drop yhat_9010

*predictedoutcome90_9010: Policy ONLY responsiven to the 90th percentile
*Note: using `coef_109010' instead of `coefntdgpt' because not dividing beta by 2 since policy only responds to one income group in DGP.
quietly: gen yhat_9010 = exp(`constant_109010' + (`coef_109010' * logodds_high) + (`coefgap_109010' * dif_90_10)) / (1 + exp(`constant_109010' + (`coef_109010' * logodds_high) + (`coefgap_109010' * dif_90_10))) if dif_90_10!=.
quietly:     gen predictedoutcome90_9010 = 0 if dif_90_10!=.
quietly:     replace predictedoutcome90_9010 = 1 if runiform() < yhat_9010 & dif_90_10!=.
drop yhat_9010

*predictedoutcome10_9010: Policy ONLY responsiven to the 10th percentile
*Note: using `coef_109010' instead of `coefntdgpt' because not dividing beta by 2 since policy only responds to one income group in DGP.
quietly: gen yhat_9010 = exp(`constant_109010' + (`coef_109010' * logodds_low) + (`coefgap_109010' * dif_90_10)) / (1 + exp(`constant_109010' + (`coef_109010' * logodds_low) + (`coefgap_109010' * dif_90_10))) if dif_90_10!=.
quietly:     gen predictedoutcome10_9010 = 0 if dif_90_10!=.
quietly:     replace predictedoutcome10_9010 = 1 if runiform() < yhat_9010 & dif_90_10!=.
drop yhat_9010

end
***************************************************************************************
***************************************************************************************



**************************************************
*Set seed 
*https://www.random.org/
set seed 45353807

//save to temporary file  
tempfile sim_difNOintercept
tempname temp1
postfile `temp1' coef290T10ngi const290T10ngi ///
			coef210T10ngi const210T10ngi ///
				sample_size  cor9010tab1 ///
			gibothinftencoef gibothinftengicoef gibothinftenconscoef gibothinftenn ///
			gibothinfntycoef gibothinfntygicoef gibothinfntyconscoef gibothinfntyn ///
			gintyinftencoef gintyinftengicoef gintyinftenconscoef gintyinftenn ///
			gintyinfntycoef gintyinfntygicoef gintyinfntyconscoef gintyinfntyn ///
			giteninftencoef giteninftengicoef giteninftenconscoef giteninftenn ///
			giteninfntycoef giteninfntygicoef giteninfntyconscoef giteninfntyn ///
				cem10ngi cem10ngicons N_cem10ngi ///
				cem90ngi cem90ngicons N_cem90ngi ///
				cem10ngiF cem10ngiconsF N_cem10ngiF ///
				cem90ngiT cem90ngiconsT N_cem90ngiT ///
				cem10ngiT cem10ngiconsT N_cem10ngiT ///
				cem90ngiF cem90ngiconsF N_cem90ngiF ///
				using "`sim_difNOintercept'", replace

//1000 runs 
forval i = 1/1000 {
*run program that simulates series
sim_90_10

****************************************************
**Generate Figures A-2 and A-3 to Confirm DGP reproduces patterns in observed data
**DGP includes group inercept
**Regression intentionally does not include group intercept
**Goal is to see if using Gilen's model recovers the observed values in the data
**when the DGP includes a group interecept.
****************************************************
*90th Percentile (policy resposive to both 10th and 90th)
quietly: logit predictedoutcome9010_9010 logodds_high  if abs(high - low)>.1
mat define b = get(_b)
gen yhat90=exp(b[1,2]+(b[1,1]*logodds_high))/(1+exp(b[1,2]+(b[1,1]*logodds_high))) if e(sample)
mat drop b

quietly: logit predictedoutcome9010_9010 logodds_high  if (high - low)>.1
mat define b = get(_b)
gen yhat90_greater=exp(b[1,2]+(b[1,1]*logodds_high))/(1+exp(b[1,2]+(b[1,1]*logodds_high))) if e(sample)
mat drop b

quietly: logit predictedoutcome9010_9010 logodds_high  if (low - high)>.1
mat define b = get(_b)
gen yhat90_lesser=exp(b[1,2]+(b[1,1]*logodds_high))/(1+exp(b[1,2]+(b[1,1]*logodds_high))) if e(sample)
mat drop b

quietly: gen high_perc = 100*high

*Save last 12 figures to use in Appendix figure
if `i' >= 989 {
quietly: twoway scatter yhat90 high_perc, msize(small) mcolor(black) ///
|| scatter yhat90_greater high_perc, mcolor(orange%80) ///
|| scatter yhat90_lesser high_perc, mcolor(navy%90) ///
lcolor(black) ylabel(, nogrid notick) xlabel(, nogrid notick) xtitle("Percent Favoring Change (Simulated 90th Income Percentile)") ///
ytitle("Predicted Probability of Policy Change") graphregion(fcolor(white) lcolor(white)) ///
plotregion(lcolor(black)) xlabel(10[10]90) xscale(r(5 95)) yscale(r(0 1)) ylabel(0(.2)1) ///
legend(off) ///
text(.085 36 "Simulated low income more" "supportive than simulated high income" "(10th - 90th)>10%") ///
text(.54 57 "Simulated high income" "more supportive" "than  simulated low income" "(90th - 10th)>10%") ///
text(.33 48 "{it:Combined}", color(black)) ///
ysize(3) xsize(4)
`FIG'graph export AppendixFigures\predictedoutcome90_9010_`i'.pdf, as(pdf) replace
}
drop yhat90 yhat90_greater yhat90_lesser


*10th Percentile (policy resposive to both 10th and 90th)
quietly: logit predictedoutcome9010_9010 logodds_low  if abs(high - low)>.1
mat define b = get(_b)
gen yhat90=exp(b[1,2]+(b[1,1]*logodds_low))/(1+exp(b[1,2]+(b[1,1]*logodds_low))) if e(sample)
mat drop b

quietly: logit predictedoutcome9010_9010 logodds_low  if (high - low)>.1
mat define b = get(_b)
gen yhat90_greater=exp(b[1,2]+(b[1,1]*logodds_low))/(1+exp(b[1,2]+(b[1,1]*logodds_low))) if e(sample)
mat drop b

quietly: logit predictedoutcome9010_9010 logodds_low  if (low - high)>.1
mat define b = get(_b)
gen yhat90_lesser=exp(b[1,2]+(b[1,1]*logodds_low))/(1+exp(b[1,2]+(b[1,1]*logodds_low))) if e(sample)
mat drop b

quietly: gen low_perc = 100*low

*Save last 12 figures to use in Appendix figure
if `i' >= 989 {
quietly: twoway scatter yhat90 low_perc, msize(small) mcolor(black) ///
|| scatter yhat90_greater low_perc, mcolor(orange%80) ///
|| scatter yhat90_lesser low_perc, mcolor(navy%90) ///
lcolor(black) ylabel(, nogrid notick) xlabel(, nogrid notick) xtitle("Percent Favoring Change (Simulated 10th Income Percentile)") ///
ytitle("Predicted Probability of Policy Change") graphregion(fcolor(white) lcolor(white)) ///
plotregion(lcolor(black)) xlabel(10[10]90) xscale(r(5 95)) yscale(r(0 1)) ylabel(0(.2)1) ///
legend(off) ///
text(.085 36 "Simulated low income more" "supportive than simulated high income" "(10th - 90th)>10%") ///
text(.54 57 "Simulated high income" "more supportive" "than  simulated low income" "(90th - 10th)>10%") ///
text(.33 48 "{it:Combined}", color(black)) ///
ysize(3) xsize(4)
`FIG'graph export AppendixFigures\predictedoutcome10_9010_`i'.pdf, as(pdf) replace
}
****************************************************
****************************************************

******************************************************************
*Group intercept excluded even though in DGP
*This is the hypothesis being tested: Observed data
*show group intercept exist. Expectation is that if this
*is the case, when analyzing preference gaps, results will
*not match DGP.
***************************
*90 vs 10
****************************************************************************


****************************************************************************
*Analysis for Table 1: Columns 2 and 4
****************************************************************************
*DGP = Responsiveness to 90th & responsivenss to 10th
quietly: logit predictedoutcome9010_9010 logodds_high  if abs(high - low)>.1
*ngi = no group intercept
mat beta = e(b)
scalar coef290T10ngi = beta[1,1]
scalar const290T10ngi = beta[1,2]

quietly: logit predictedoutcome9010_9010 logodds_low  if abs(high - low)>.1
mat beta = e(b)
scalar coef210T10ngi = beta[1,1]
scalar const210T10ngi = beta[1,2]
//save sample size 
scalar sample_size = e(N)

//Correlation between simulated preferences when preference gap exists
quietly: cor high low if abs(high - low)>.1
scalar cor9010tab1 = r(rho)


****************************************************************************
*Analysis for Appendix Table A-5
*Test whether group intercept recovers correct DGP
****************************************************************************
**10th vs 90: Both Influence DGP (predictedoutcome9010_9010)
quietly: logit predictedoutcome9010_9010 logodds_low dif_90_10
mat def coef = e(b)
scalar gibothinftencoef = coef[1,1]
scalar gibothinftengicoef = coef[1,2]
scalar gibothinftenconscoef = coef[1,2]
scalar gibothinftenn = e(N)

quietly: logit predictedoutcome9010_9010 logodds_high dif_90_10
mat def coef = e(b)
scalar gibothinfntycoef = coef[1,1]
scalar gibothinfntygicoef = coef[1,2]
scalar gibothinfntyconscoef = coef[1,2]
scalar gibothinfntyn = e(N)

**10th Does NOT influence DGP, 90th does (predictedoutcome90_9010)
quietly: logit predictedoutcome90_9010 logodds_low dif_90_10
mat def coef = e(b)
scalar gintyinftencoef = coef[1,1]
scalar gintyinftengicoef = coef[1,2]
scalar gintyinftenconscoef = coef[1,2]
scalar gintyinftenn = e(N)

quietly: logit predictedoutcome90_9010 logodds_high dif_90_10
mat def coef = e(b)
scalar gintyinfntycoef = coef[1,1]
scalar gintyinfntygicoef = coef[1,2]
scalar gintyinfntyconscoef = coef[1,2]
scalar gintyinfntyn = e(N)

**10th Does influence DGP, 90th does NOT (predictedoutcome10_9010)
quietly: logit predictedoutcome10_9010 logodds_low dif_90_10
mat def coef = e(b)
scalar giteninftencoef = coef[1,1]
scalar giteninftengicoef = coef[1,2]
scalar giteninftenconscoef = coef[1,2]
scalar giteninftenn = e(N)

quietly: logit predictedoutcome10_9010 logodds_high dif_90_10
mat def coef = e(b)
scalar giteninfntycoef = coef[1,1]
scalar giteninfntygicoef = coef[1,2]
scalar giteninfntyconscoef = coef[1,2]
scalar giteninfntyn = e(N)
****************************************************************************
****************************************************************************


****************************************************************************
*Analysis for Table 3: 
*Test whether CEM recovers correct DGP
****************************************************************************
quietly: drop if dif_90_10==.

*10th vs 90: Both Influence DGP (predictedoutcome9010_9010)
*10th income percentile
quietly: cem low if predictedoutcome9010_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome9010_9010 logodds_low if cem_matched==1

mat def coef = e(b)
scalar cem10ngi = coef[1,1]
scalar cem10ngicons = coef[1,2]
scalar N_cem10ngi = e(N)

*90th income percentile
quietly: cem high if predictedoutcome9010_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome9010_9010 logodds_high if cem_matched==1

mat def coef = e(b)
scalar cem90ngi = coef[1,1]
scalar cem90ngicons = coef[1,2]
scalar N_cem90ngi = e(N)

*10th Does NOT influence DGP, 90th does (predictedoutcome90_9010)
*10th income percentile
quietly: cem low if predictedoutcome90_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome90_9010 logodds_low if cem_matched==1

mat def coef = e(b)
scalar cem10ngiF = coef[1,1]
scalar cem10ngiconsF = coef[1,2]
scalar N_cem10ngiF = e(N)

*90th income percentile
quietly: cem high if predictedoutcome90_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome90_9010 logodds_high if cem_matched==1
mat def coef = e(b)
scalar cem90ngiT = coef[1,1]
scalar cem90ngiconsT = coef[1,2]
scalar N_cem90ngiT = e(N)

*10th Does influence DGP, 90th does NOT (predictedoutcome10_9010)
*10th income percentile
quietly: cem low if predictedoutcome90_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome90_9010 logodds_low if cem_matched==1

mat def coef = e(b)
scalar cem10ngiT = coef[1,1]
scalar cem10ngiconsT = coef[1,2]
scalar N_cem10ngiT = e(N)

*90th income percentile
quietly: cem high if predictedoutcome90_9010!=., tr(dif_90_10) k2k
quietly: logit predictedoutcome90_9010 logodds_high if cem_matched==1
mat def coef = e(b)
scalar cem90ngiF = coef[1,1]
scalar cem90ngiconsF = coef[1,2]
scalar N_cem90ngiF = e(N)


****************************************************************************
****************************************************************************


*****************************************************************************  
//post results to dataset 
post `temp1' (coef290T10ngi) (const290T10ngi) ///
			(coef210T10ngi) (const210T10ngi) ///
			(sample_size) (cor9010tab1) ///
			(gibothinftencoef) (gibothinftengicoef) (gibothinftenconscoef) (gibothinftenn) ///
			(gibothinfntycoef) (gibothinfntygicoef) (gibothinfntyconscoef) (gibothinfntyn) ///
			(gintyinftencoef) (gintyinftengicoef) (gintyinftenconscoef) (gintyinftenn) ///
			(gintyinfntycoef) (gintyinfntygicoef) (gintyinfntyconscoef) (gintyinfntyn) ///
			(giteninftencoef) (giteninftengicoef) (giteninftenconscoef) (giteninftenn) ///
			(giteninfntycoef) (giteninfntygicoef) (giteninfntyconscoef) (giteninfntyn) ///
			(cem10ngi) (cem10ngicons) (N_cem10ngi) ///
			(cem90ngi) (cem90ngicons) (N_cem90ngi) ///
			(cem10ngiF) (cem10ngiconsF) (N_cem10ngiF) ///
			(cem90ngiT) (cem90ngiconsT) (N_cem90ngiT) ///
			(cem10ngiT) (cem10ngiconsT) (N_cem10ngiT) ///
			(cem90ngiF) (cem90ngiconsF) (N_cem90ngiF) 
}

mat l betaDGP

postclose `temp1' //write results to temporary dataset

use "`sim_difNOintercept'", clear //load new dataset with simulated data  

*******************************


local varlist coef290T10ngi const290T10ngi coef210T10ngi const210T10ngi ///
sample_size cor9010tab1 ///
gibothinftencoef gibothinftengicoef gibothinftenconscoef gibothinftenn ///
gibothinfntycoef gibothinfntygicoef gibothinfntyconscoef gibothinfntyn ///
gintyinftencoef gintyinftengicoef gintyinftenconscoef gintyinftenn ///
gintyinfntycoef gintyinfntygicoef gintyinfntyconscoef gintyinfntyn ///
giteninftencoef giteninftengicoef giteninftenconscoef giteninftenn ///
giteninfntycoef giteninfntygicoef giteninfntyconscoef giteninfntyn ///
cem10ngi cem10ngicons N_cem10ngi ///
cem90ngi cem90ngicons N_cem90ngi ///
cem10ngiF cem10ngiconsF N_cem10ngiF cem90ngiT cem90ngiconsT N_cem90ngiT  ///
cem10ngiT cem10ngiconsT N_cem10ngiT cem90ngiF cem90ngiconsF N_cem90ngiF  //local list of all variables 

foreach var of local varlist { //generate locals of means and 95% bounds
 	_pctile `var', nq(1000)
	local lb_`var' = r(r25) 
	local m_`var' = r(r500) 
    local ub_`var' = r(r975) 	
}

//export results 


*********************
*Table 1
*********************
*Col 2
di `m_coef210T10ngi',`lb_coef210T10ngi',`ub_coef210T10ngi'
di `m_const210T10ngi',`lb_coef210T10ngi',`ub_const210T10ngi'
di `m_sample_size'
*Col 4
di `m_coef290T10ngi',`lb_coef290T10ngi',`ub_coef290T10ngi'
di `m_const290T10ngi',`lb_const290T10ngi',`ub_const290T10ngi'
di `m_sample_size'
*********************
*nt = 90 True
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ntngi) result(`m_coef290T10ngi') append // 90th percentile beta (responsiveness)
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbntngi) result(`lb_coef290T10ngi') append //lower bound
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubntngi) result(`ub_coef290T10ngi') append //upper bound

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cngi) result(`m_const290T10ngi') append // constant
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcngi) result(`lb_const290T10ngi') append //lower bound
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcngi) result(`ub_const290T10ngi') append //upper bound
*tt = 10 true
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ttngi) result(`m_coef210T10ngi') append //10th percentile beta (responsiveness)
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbttngi) result(`lb_coef210T10ngi') append // lower bound
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubttngi) result(`ub_coef210T10ngi') append // upper bound

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ctngi) result(`m_const210T10ngi') append // constant
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbctngi) result(`lb_const210T10ngi') append //lower bound
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubctngi) result(`ub_const210T10ngi') append //upper bound

**sample size mean 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(nngi) result(`m_sample_size') round(0) append // average sample size across simulations
*Correlation b/t 90th and 10th income percentile when preference gap exists 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(corpreftabone) result(`m_cor9010tab1') round(2) append //

********************
*Table 3
********************
*tt = 10 true CEM 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngi) result(`m_cem10ngi') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngi) result(`lb_cem10ngi') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngi) result(`ub_cem10ngi') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngicons) result(`m_cem10ngicons') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngicons) result(`lb_cem10ngicons') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngicons) result(`ub_cem10ngicons') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemtngi) result(`m_N_cem10ngi') round(0) append
*nt = 90 True CEM
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngi) result(`m_cem90ngi') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngi) result(`lb_cem90ngi') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngi) result(`ub_cem90ngi') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngicons) result(`m_cem90ngicons') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngicons) result(`lb_cem90ngicons') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngicons) result(`ub_cem90ngicons') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemnngi) result(`m_N_cem90ngi') round(0) append

* = 10 false CEM 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngiF) result(`m_cem10ngiF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngiF) result(`lb_cem10ngiF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngiF) result(`ub_cem10ngiF') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngiconsF) result(`m_cem10ngiconsF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngiconsF) result(`lb_cem10ngiconsF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngiconsF) result(`ub_cem10ngiconsF') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemtngiF) result(`m_N_cem10ngiF') round(0) append
*nt = 90 (only) True CEM
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngiT) result(`m_cem90ngiT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngiT) result(`lb_cem90ngiT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngiT) result(`ub_cem90ngiT') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngiconsT) result(`m_cem90ngiconsT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngiconsT) result(`lb_cem90ngiconsT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngiconsT) result(`ub_cem90ngiconsT') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemnngiT) result(`m_N_cem90ngiT') round(0) append

* = 10 (only) True CEM 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngiT) result(`m_cem10ngiT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngiT) result(`lb_cem10ngiT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngiT) result(`ub_cem10ngiT') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemtngiconsT) result(`m_cem10ngiconsT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemtngiconsT) result(`lb_cem10ngiconsT') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemtngiconsT) result(`ub_cem10ngiconsT') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemtngiTt) result(`m_N_cem10ngiT') round(0) append
*nt = 90 False CEM
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngiF) result(`m_cem90ngiF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngiF) result(`lb_cem90ngiF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngiF) result(`ub_cem90ngiF') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(cemnngiconsF) result(`m_cem90ngiconsF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbcemnngiconsF) result(`lb_cem90ngiconsF') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubcemnngiconsF) result(`ub_cem90ngiconsF') append

`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(mNcemnngiFn) result(`m_N_cem90ngiF') round(0) append

**********************
*Appendix Table A-5: Does including group intercept in model recover correct dgp?
**********************
*both 10th and 90th influence policy outcomes
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinftencoef) result(`m_gibothinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinftengicoef) result(`m_gibothinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinftenconscoef) result(`m_gibothinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinftenn) result(`m_gibothinftenn') round(0) append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinfntycoef) result(`m_gibothinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinfntygicoef) result(`m_gibothinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinfntyconscoef) result(`m_gibothinfntyconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gibothinfntyn) result(`m_gibothinfntyn') round(0) append
*lb
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinftencoef) result(`lb_gibothinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinftengicoef) result(`lb_gibothinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinftenconscoef) result(`lb_gibothinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinfntycoef) result(`lb_gibothinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinfntygicoef) result(`lb_gibothinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgibothinfntyconscoef) result(`lb_gibothinfntyconscoef') append
*ub
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinftencoef) result(`ub_gibothinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinftengicoef) result(`ub_gibothinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinftenconscoef) result(`ub_gibothinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinfntycoef) result(`ub_gibothinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinfntygicoef) result(`ub_gibothinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgibothinfntyconscoef) result(`ub_gibothinfntyconscoef') append
*Only 90th influence policy outomces
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinftencoef) result(`m_gintyinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinftengicoef) result(`m_gintyinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinftenconscoef) result(`m_gintyinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinftenn) result(`m_gintyinftenn') round(0) append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinfntycoef) result(`m_gintyinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinfntygicoef) result(`m_gintyinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinfntyconscoef) result(`m_gintyinfntyconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(gintyinfntyn) result(`m_gintyinfntyn') round(0) append
*lb
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinftencoef) result(`lb_gintyinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinftengicoef) result(`lb_gintyinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinftenconscoef) result(`lb_gintyinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinfntycoef) result(`lb_gintyinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinfntygicoef) result(`lb_gintyinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgintyinfntyconscoef) result(`lb_gintyinfntyconscoef') append
*ub
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinftencoef) result(`ub_gintyinftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinftengicoef) result(`ub_gintyinftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinftenconscoef) result(`ub_gintyinftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinfntycoef) result(`ub_gintyinfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinfntygicoef) result(`ub_gintyinfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgintyinfntyconscoef) result(`ub_gintyinfntyconscoef') append
*Only 10th influences policy outcomes
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninftencoef) result(`m_giteninftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninftengicoef) result(`m_giteninftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninftenconscoef) result(`m_giteninftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninftenn) result(`m_giteninftenn') round(0) append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninfntycoef) result(`m_giteninfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninfntygicoef) result(`m_giteninfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninfntyconscoef) result(`m_giteninfntyconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(giteninfntyn) result(`m_giteninfntyn') round(0) append
*lb
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninftencoef) result(`lb_giteninftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninftengicoef) result(`lb_giteninftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninftenconscoef) result(`lb_giteninftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninfntycoef) result(`lb_giteninfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninfntygicoef) result(`lb_giteninfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(lbgiteninfntyconscoef) result(`lb_giteninfntyconscoef') append
*ub 
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninftencoef) result(`ub_giteninftencoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninftengicoef) result(`ub_giteninftengicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninftenconscoef) result(`ub_giteninftenconscoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninfntycoef) result(`ub_giteninfntycoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninfntygicoef) result(`ub_giteninfntygicoef') append
`TEX'texresults2 using TxtFiles_NumericalResults\results_simngint.txt, texmacro(ubgiteninfntyconscoef) result(`ub_giteninfntyconscoef') append






	